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ABSTRACT 


Transient analysis of large structural systems is a computationally demanding 
process, which in the past has prevented dynamic redesign and optimization. Large 
structures, such as buildings or ships subjected to random base motions, use isolators to 
minimize strain energies, which may cause damage or structural failure. 

This research focuses on the optimization of isolator parameters in order for 
structural systems to withstand potentially catastrophic transient vibrations. Many non¬ 
linear hysteretic, viscoelastic, and sliding friction isolators were numerically modeled 
using the scientific programming language, MATLAB. The existing programs used to 
solve the Voltera integral formulation for Transient Structural Synthesis (TSS) and the 
Recursive Block-by-Block (RBB) algorithm were investigated and enhanced to yield 
greater accuracy and increased computational speed. The final product is a user-friendly 
Decision Support System (DSS) for use with both civil and military applications. Based 
on different types of base motions and the inherent dynamics of the structure, this DSS is 
capable of optimizing isolator parameters to meet a user specific objective. 
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EXECUTIVE SUMMARY 


On average, there are about 500,000 deteetable earthquakes in the world per year. 
Only about 100,000 of these ean be felt and about 100 of them have caused moderate to 
severe damage. Earthquakes strike without warning and can cause wide-spread damage 
in populated areas, particularly in third world countries where building codes are poorly 
enforced or non-existent. 

Shear failure is the predominant failure mechanism of most structural systems 
experiencing a seismic event. Passive isolators have proved successful in absorbing the 
shear energy delivered by strong motion earthquakes. However, selection of the correct 
type of isolator that will respond best to the inherent dynamic characteristics of a specific 
structural system poses the greatest challenge. Traditional methods of transient analysis, 
i.e. Finite Element Analysis (FEA), for complex structural systems are computationally 
demanding and make it impracticable for use with existing optimization techniques. 

This research focuses on the optimization of isolator parameters in order for 
structural systems to withstand potentially catastrophic transient vibrations. Many non¬ 
linear hysteretic, viscoelastic, and sliding friction isolators were numerically modeled 
using the scientific programming language, MATLAB. The existing programs used to 
solve the Voltera integral formulation for Transient Structural Synthesis (TSS) and the 
Recursive Block-by-Block (RBB) algorithm were investigated and enhanced to yield 
greater accuracy and increased computational speed. The final product is a user-friendly 
Decision Support System (DSS) for use with both civil and military applications. Based 
on different types of base motions and the inherent dynamics of the structure, this DSS is 
capable of optimizing isolator parameters to meet a user specific objective. 
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I. INTRODUCTION 


The U.S. Geological Survey estimates an average of 500,000 detectable 
earthquakes in the world per year [Ref 1]. Only about 100,000 of these can be felt and 
about 100 of them have caused moderate to severe damage. Earthquakes strike without 
warning and can cause wide-spread damage in populated areas, particularly in third world 
countries where building codes are poorly enforced or non-existent. 

Over the past century, scientists and engineers have investigated methods to 
isolate large structural systems from strong-base motion. The large lateral displacements 
associated with strong motion flood structural systems with high-magnitude shear forces. 
Shear failure is the predominant failure mechanism of most structural system undergoing 
an earthquake event. Passive isolators have been proven successful in absorbing the 
shear energy delivered by earthquakes [Ref. 2]. 

Even though passive structural isolators present solutions to the problem of 
partially decoupling structural systems from ground motion, selection of the correct type 
of isolators that respond best to the inherent dynamic characteristics of distinct structural 
systems poses the greatest challenge. This type of problem is clearly one of optimization 
- choose the best set of parameters to minimize or maximize an objective. Traditional 
methods of transient analysis, i.e. Einite Element Analysis (EEA), for complex structural 
systems are computationally demanding. They are impracticable for use with existing 
optimization techniques due to the computational demand placed by nonlinear isolators. 

Transient Structural Synthesis (TSS) is a newer transient analysis method that 
uses normalized mode shapes, natural frequencies, and an integral equation formulation 
to solve the response of complex structural systems with localized nonlinearities [Ref. 3]. 
This time-domain synthesis method uses the Voltera Integral equation, which relies on 
the exact solution to the governing equations of motion. This method greatly decreases 
processing speed by allowing implicit model reduction and exploiting the contractive 
nature of the Voltera Integral equation. 

Eurther improvements in processing speed can be achieved with the use of a 
Recursive Block-by-Block algorithm [Ref 14]. This method divides the simulation time 


1 



and the impulse response function solutions into equally-sized blocks. The blocks are 
individually solved by convolution and retained. The solutions to the proceeding blocks 
are dependent on the retained solutions of previously processed blocks. This method 
does not degrade the solution’s accuracy when compared with non-blocked solutions 

These methods present excellent savings in computation speed, but they were not 
without fault. The existing MATLAB programs that executed the TSS and RBB 
algorithms had programming redundancies that decreased processing time. They also 
yielded less accurate results since they employed MATLAB’s conv command, which 
uses a rectangular integration rule. The non-linear isolator models developed from a 
previous thesis research did not execute correctly with the RBB algorithm, which 
introduced errors into the final synthesized solutions. 

Presented here is a thorough investigation of the existing MATLAB programs, 
which led to optimized program logic, corrected non-linear isolator models, and a method 
that transforms MATLAB’s rectangular integration rule convolution into any integration 
rule with higher-order error. Also presented are numerical models for sliding-friction 
isolators, and a user friendly Decision Support System (DSS) tool for transient response 
and analysis. This program gives a user the ability to quickly change isolator parameters 
for TSS / RBB, change types of base motions, and optimize the isolator parameters to 
meet user defined objectives. 
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II. BACKGROUND 


A. EARTHQUAKES 

Earthquakes are one of the most terrifying natural phenomena known. They 
indiseriminately strike without warning, are felt thousands of miles away, and ean cause 
widespread damage. The largest earthquake to strike in the twentieth century was in 
1960 off the coast of Chile - 9.5 on the Richter scale. It ruptured an area 525 miles long 
and 80 miles wide, which devastated numerous Chilean cities, caused hundred of large 
landslides, and may have aggravated the Puyehue Volcano that erupted two days later 
[Ref 4]. This earthquake created a deadly tsunami of more than 30 feet that completely 
eliminated entire villages along the Chilean coast. The tsunami continued across the 
Pacific Ocean where it killed 61 people in Hawaii and hundreds more in Japan, more than 
8000 miles from the earthquake source. 

The causes of earthquakes are explained with the theory of plate tectonics. The 
basic theory is that the lithosphere, the surface layer of the earth, is made up of many 
plates that slide over the asthenoshpere, the layer that behaves as a lubricant. Along the 
boundaries of the plates they will move apart, push together, or slide against each other. 
When two plates move apart from each other, molten lava rises from beneath the 
lithosphere. The lava fills the gap created from the separating plates and hardens. This 
divergent plate boundary is commonly seen at the bottom of an ocean. 

Along convergent plate boundaries, one of the plates is usually pushed under the 
other. The subducted plate will sink into the lower mantel layers and melt. When neither 
of the plates yields to subduction, the colliding plates form mountain ranges. The peak of 
Mt. Everest, part of the Himalaya mountain range, increases an average of three inches 
per year due to the pressure created from the Asian Plate and the Indian Plate [Ref. 5]. 
Plates that slide against each other do so along transform boundaries and can create 
tension. The sudden decrease in tension or compression between moving plates releases 
large amounts of shear energy that propagate through the Earth’s crust. Earthquakes are 
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more common along fault lines - areas where plates meet and move in different 
directions. The plate-boundary diagram presented in Figure 1 further explains these 
phenomena. 




CONVERQENT TRANSFORM DIVERGENT CONVERGENT 

PLATE BOUNDARY BOUNDARY PLATE BOUNDARY PLATE BOUNDARY 


CONTINENTAL RlFTZONE 
(YOUNG PLATE BOUNDARY) 



f/nic s^readino 

RiOCC 


LITHOSPHERE 



ASTHENOSPHERE 


SUCrvCTlNG 

P.ATE 


Figure 1. Plate-Boundary Diagram (From Ref 1) 


1, Fault Lines 

Geologists have identified four types of fault systems, normal, reverse, thrust, and 
strike-slip faults. They are characterized by orientation of the fault plane, the break in the 
rock, and the resulting movement of the two blocks of rock. The normal fault has a 
nearly vertical fault plane with a hanging wall, the rock above the plane, and a foot wall, 
the rock below the plane. The hanging wall pushes down against the foot wall, which in 
turn pushes up against the hanging wall. This separation of the Earth’s crust is a result of 
the pull from a divergent plate boundary. The reverse fault also has a nearly vertical fault 
plane but the compressed rock causes the hanging wall to push up and the footwall to 
push down. A thrust fault behaves similar to the reverse fault except that the fault plane 
is nearly horizontal. As a result, the hanging wall is pushed on top of the foot wall - a 
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common occurrence in a eonvergent plate boundary. The strike-slip fault commonly 
occurs in transform plate boundaries where the plates move in opposing horizontal 
direetions. 


Normal Fault 



Figure 2. Normal Fault (From Ref. 6) 


Reverse Fault 



Figure 3. Reverse Fault (From Ref. 6) 
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Figure 4. Strike-Slip Fault (From Ref. 6) 


The movement of the plates along fault lines ereates frietion. When the moving 
rooks beoome looked, the restricted motion causes the build up of potential energy - the 
stress in the rock increases. When the rock finally yields due to the enormous stresses, 
the sudden release of kinetic energy propagates through the Earth’s lithosphere. 

2, Seismic Waves 

The seismic energy wave propagates through the Earth in the same manner as 
disturbance in a body of water. Seismic waves are categorized into body waves and 
surface waves. Body waves are those that travel through the inner part of the Earth and 
cause primary waves (p-waves) and secondary waves (s-waves) [Ref 6]. P-waves are 
compression waves that travel through to the Earth’s surface the fastest and can penetrate 
solids, liquids, and gasses. S-waves are shear waves that travel perpendicular to the 
direction of the body waves. They only propagate through solid material and are thus 
stopped by the Earth’s liquid asthenosphere. 
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Focus 


Figure 5. Body Wave Propagation (From Ref. 6) 

Surface waves are caused by the body waves that reach the surface. They are also 
known are long waves (1-waves) due to their low frequency. L-waves produce the most 
intense ground displacement and are responsible for the most damage associated with 
earthquakes. 

B. BASE MOTION ISOLATION 

Base motion isolation has become an increasing requirement for structural designs 
in earthquake prone areas. The two types of isolation that have gained acceptance are 
elastomeric bearings and sliding systems. Elastomeric bearings are made of either 
natural rubber or neoprene and serve the purpose of decoupling a structural system from 
the horizontal component of strong-ground motion. The low horizontal stiffness of the 
elastomer gives the structural system a natural frequency that is much lower than its 
fixed-base frequency and lower than the larger ground motion frequencies [Ref 2]. 
Since the first non-rigid body mode of a structural system relates to structural 
deformations, the system can be viewed as a rigid body. 

Sliding systems limit the transfer of shear across the isolation interface. The 
sliding interface is typically made of two dissimilar metals with a low friction coefficient 
aided by a restoring device. The use of restoring devices with sliding systems helps 
reduce the diameter of the bearing plates and the support system. 
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Non-linear base isolation is of speeial interest sinee, unlike linear isolators that 
defleet earthquake energy; non-linear isolators dissipate the energy through plastic 
deformation of the isolator materials or with the use of viscous or coulomb damping. 

1, Elastomeric and Hysteretic Isolators 

Elastomeric isolators usually include high-damping rubber bearings (HDR), low- 
damping rubber bearings (RB) or low-damping rubber bearings with a lead core (LRB). 
These isolators consist of thin layers or natural rubber that are vulcanized and bonded to 
steel plates. Natural rubber exhibits both hysteretic and viscoelastic behavior. In 
particular, low-damping natural rubber exhibits an almost linear elastic and linear viscous 
behavior. HDR and RB isolators are practically linear since their dynamic shear modulus 
remains almost constant over a wide range of temperature and frequencies [Ref. 2]. Their 
dynamic damping coefficients remain comparatively constant at frequencies lower than 
50 Hz - a desirable characteristic since strong motion frequencies are much lower. 

For isolators with different types of rubber, the stiffness characteristics will 
change as the shear loading changes. For small displacements, the isolator stiffness is 
large - this is desirable to minimize the effects of wind loading. When the isolator 
experiences large shear strains, the stiffness is very low and effectively deflects seismic 
energy. 



Figure 6. Idealized Hysteretic Force-Displacement Relation of a FRB (From Ref. 7) 
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Bilinear Force - Displacennent Loop 



Figure 7. Bilinear Hysteretie Numerieal Model of a LBR 

A LRB is eategorized as a hysteretie isolator beeause of its lead eore. A 
hysteretie material is one that will remain permanently deformed when the disturbing 
force is removed. Viscoelastic and elastomeric materials do not exhibit hysteretie 
behavior since they return to their original position once the disturbing force is removed. 
LRB isolators are typically constructed of low-damping natural rubber with a press fitted 
lead core. During strong lateral motion, the lead core yields at low stresses, deforms in 
nearly pure shear, and produces stable hysteretie behavior. Repeated yielding will not 
cause fatigue failure since lead self-anneals at standard temperature [Ref. 7]. The 
behavior of LRB isolators can be represented by a bilinear hysteretie numerical model, 
such as the one presented in Figure 7. Equation (1) relates the characteristic strength, Q, 
of an LRB isolator to the lead core area, Ap, and the shear yield stress of lead, Qyi. 

Q = ApCr^, ( 1 ) 

Even though a bilinear numerical model will yield acceptable results for non¬ 
linear isolators, the Wen model, Eigure 8, will more closely represent the hysteretie 
curve. 
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Wen Hysteresis Loop 



2, Sliding Isolators with Restoring Force Mechanisms 

Sliding bearings will limit the transmission of lateral forces to the isolated 
structure, but require large sliding surfaces and usually result in permanent offset 
displacements. To reduce the diameter of the sliding surface and prevent unwanted 
displacements, a restoring force mechanism is typically used. The three types of sliding 
isolators are shown in Figure 9. 


Flat sliding 
surface 
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Figure 9. Idealized Force Displacement Loops of Sliding Bearings (From Ref 7) 
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The lateral force developed from a sliding bearing can be calculated with 
Equation (2). The effects due to vertical ground acceleration and overturning moments 
on the normal load used to calculate the lateral force is calculated with Equation (3) [Ref 
V]. 

F = ^U + M^.Nsgn{U) ( 2 ) 

K 

where 

U = Displacement 
Udot = Sliding velocity 

R = Radius of curvature of the sliding surface 
Ps = Coefficient of sliding friction 
N = Normal load on the bearing 


N = W 


^, U 

1 + ^ + ^ 
g w 


( 3 ) 


Elat sliding bearings have an infinite radius of curvature and require a restoring 
mechanism. This is not the case with spherical and conical sliding surfaces since self¬ 
centering is inherent to their design. Sliding bearings are commonly made from PTEE or 
PTEE based composites in contact with a polished stainless steel surface. 

3, Viscoelastic Isolators 

Viscoelastic isolators can be either frequency dependent or frequency 
independent. Such isolators that resemble a parallel spring-damper configuration are 
frequency dependent and poorly dissipate displacement energy. In turn, if the 
configuration resembles a spring-damper in series, such as the Maxwell model, it is 
considered frequency independent. Erequency independent isolators retain their 
characteristic force-displacement loop over a broad range of frequencies [Ref 8]. 

C. RELEVANT RESEARCH 

The final analysis / optimization program is a combination of the proper use of the 
Voltera Integral Equation, Transient Structural Synthesis (TSS), and the Recursive 


11 



Block-by-Block (RBB) algorithm. These topies are well researehed and diseussed in 
their respeetive referenees. The following sub-seetions broadly diseuss these topies in 
order for the reader to have the general knowledge need to understand the following 
ehapters. 

1. Governing Equations of Motion 

Meehanieal or struetural systems undergo foreed vibrations when external energy 
is applied. External energy ean be applied by an external / internal foree or by a ehange 
in displaeement. Forees and displaeements ean be harmonie, non-harmonie and periodie, 
non-periodie, or random in nature. Harmonie exeitation usually results in a harmonie 
system response, while suddenly applied non-periodie exeitation results in a transient 
system response. However, suddenly applied harmonie exeitation will also produee a 
transient response. 

The basie equation of motion for a foreed single degree-of-freedom (DOF) system 
is derived from Newton’s seeond law of motion. 


mx + cx + kx = F{t) +F^^{x,x,t) (4) 

The solution to this order ordinary differential equation is a homogeneous and 
partieular one. When the initial eonditions are zero, sueh as a struetural system awaiting 
a seismie event, the homogeneous solution will also be zero. What remains is the 
partieular solution to the foreed system. Finear systems ean be solved with the eorreet 
Faplaee transform, while nonlinear systems may require more advaneed numerieal 
methods. 

A eommon method for solving linear-multi DOF systems is with eigenveetor and 
eigenvalue solutions. For an undamped 2-DOF system. Figure 10, the equations of 
motion are stated here in matrix form. 



0 ^ 


4- 


J 







vO 

m^j 

1^2 J 

1 -K 

K j 

1 



( 5 ) 
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Figure 10. Forced, 2-DOF Undamped System 

The mass and stiffness matrices are needed in order to solve for the eigenvalues 
and eigenvectors. The square root of the eigenvalues yields the natural frequencies and 
the eigenvectors contain the system modes. Systems with n-DOF will have n-modes and 
n-natural frequencies. The modes and natural frequencies are needed to later calculate 
the impulse response functions that will be used to solve the particular solution using the 
Voltera integral equation. 


2. Voltera Integral Equation 

This convolution integral is generally found in most college-level vibrations 
textbooks and is based on the idea that any force is the equivalent to a superposition of 
impulses of discrete duration time, x. 

= ( 6 ) 


Here, h are the impulse response functions (IRF) that were calculated with the 
modes and natural frequencies. The equation used to calculate the IRFs for a 1-DOF 
system is: 


h{t) = 




sin(y„ 


y 


(V) 


where 

cOn = Natural frequency 
C, = Damping ratio 

In order to solve the convolution integral numerically, an integration rule is 
applied in matrix form. The following matrix equation was developed with the trapezoid 
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integration rule and solves for the partieular solution of a system with zeros for the initial 
eonditions. 




(OAt) 
-h (lAO 


1 


h OAt 


(mAt) h((m-l)At) ••• (OAt) 


T^ISO 

^0 


FI: 


( 8 ) 


The forcing vector represents the isolator force due to the base motion 
excitation. Reference [8] presents a detailed discussion of the Voltera Integral Equation 
in modal space. 

3, Transient Structural Synthesis 

Transient Structural Synthesis (TSS) is a topic that is thoroughly discussed in 
[Ref 3] and broadly presented here since it is the core of the program developed for this 
thesis. In a time domain, TSS is used to calculate the dynamic response of coupled 
substructures and allows for model reduction to only the coordinates of interest. System 
modifications can quickly be calculated with pre-modification transient response data. 
The object is to solve for the nonlinear isolator forces, which are retained to later solve 
for the response of other coordinates of interest. This is an iterative process that requires 
an initial guess for f®", the vector of isolator forces - usually a vector of ones. Chapter 
III, Section C presents a method of improving the initial guess for accuracy by modifying 
f*” with a vector of integration rule weightings. The following flowchart explains this 
iterative method. 
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Figure 11. Iterative method used to solve for the isolator forces 


The converged solution for the vector of isolator forces is used to calculate the 
dynamic response of other coordinates of interest. The response of such a point is the 
linear combination of its IRF convolved with each isolator force that couples the system. 

4, Recursive Block-by-Block Algorithm 

The recursive block-by-block (RBB) algorithm [Ref 14] further reduces 
processing time by dividing the total simulation time into several blocks of equal lengths. 
This method divides the simulation time and the impulse response function solutions into 
equally sized blocks. The blocks are individually solved by convolution and retained. 
The solutions to the proceeding blocks are dependent on the retained solutions of 
previously processed blocks. This method does not degrade the solution’s accuracy when 
compared with non-blocked solutions. 

D, OPTIMIZATION - NONLINEAR PROGRAMMING 

Optimization problems can be divided into either linear programming (LP) or 
nonlinear programming (NLP) - most engineering problems are nonlinear in nature. 
NLP algorithms start at any point within the feasible solution region of the NLP problem. 
It then attempts to move from the starting point in a direction through the feasible region 
that results in an improved objective function. If the algorithm determines there is a 
better direction, it will adjust its step size and move in that direction to a better feasible 
solution. The process continues until the algorithm determines that movement in any 
direction will degrade the solution. The algorithm terminates when it has found a local 
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optimal solution. - a solution that is better than any other feasible solution in its 
immediate area within the feasible region. However, a given loeal optimal solution might 
not be the best possible solution. It is eommon praetiee with NLP problems to use 
different starting points in seareh of the best feasible solution [Ref. 10]. 



Figure 12. Multiple Solutions to a NLP Problem (From Ref 10) 


NLP problems have three eombinations of objeetives with eonstraints: linear 
objeetive with nonlinear eonstraints, nonlinear objeetive with linear eonstraints and 
nonlinear objeetive with nonlinear eonstraints. The Deeision Support System (DSS) tool 
written for this thesis uses a linear objeetive funetion of either displaeements, 
aeeelerations, or both, with nonlinear eonstraints. These eombinations are graphieally 
deseribed in Figure 13. 
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Figure 13. Different Combinations of Objeetives with Constraints (From Ref 10) 


1, One-Variable Optimization Problems 

For multivariable funetions, a one-dimensional minimization routine may be 
ealled numerous times to find a loeal minimum. Interval bounding and Golden Seetion 
interval refinement algorithms are useful in searehing for the minimum of one-variable 
funetions. 

The interval bounding funetion starts at a point on a eurve and “marehes” at a 
fixed spaee step in a deseending direetion. When the funetion eneounters an aseending 
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direction, it stops and saves the current and last two points for the interval refinement 
function. This function should be run several times at different starting points to find the 
interval that may yield an optimal solution. 

The interval refinement function uses the Golden Section algorithm to isolate the 
contained minimum point to a predetermined tolerance. This algorithm uses the Golden 
Ratio (0.618034) to refine an interval. It first determines the interval length and then 
adds a fraction of the interval length (the length multiplied by the Golden Ratio) to the 
lower bound of the interval. It subtracts the same length fraction from the upper bound. 
Depending on the magnitude of the two new points within the interval, the algorithm 
eliminates the portion of the interval that is no longer needed - Figure 14 clarifies this 
methodology [Ref 11]. 



Figure 14. Golden Section Interval Refinement (From Ref. 11) 


2. Multi-Variable Optimization Problems 

Common methods for multi-variable unconstrained optimization include Steepest 
Descent, Newton’s Method, Marquardt’s Compromise, and Modified Newton’s Methods. 
The Steepest Descent and Modified Newton methods incorporate the bounding phase and 
Golden Section algorithms. The MATLAB Optimization Toolbox used for this thesis 
automatically chooses a method; however, a user can easily change the method and 
optimization parameters by changing the optimset values within MATLAB. The 
advantages and disadvantages of each method are summarized in Table 1. 
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Method 

Advantages 

Disadvantages 

Steepest Deseent 

-Maintains a deseent property 
-Only uses the gradient 
information 

-Converges from any starting point 
-Reliable and stable 

-Requires a line seareh 
-Slows down near a minimum 
-Very slow eonvergenee 

Newton 

-Can exhibit quadratie 
eonvergenee 
-No line seareh 
-Very fast eonvergenee 

-Requires the Hessian and its inverse 
(will fail if singular) 

-Needs to solve a linear system of 
equations 

-Only eonverges near optimal 
funetion values 

-May diverge for non-quadratie 
funetions or with poor starting 
points 

Marquardt 

-Maintains a deseent property 
-No line seareh 
-Exeellent eonvergenee rate 
-Uses Newton when it gives good 
deseent results and Steepest 
Deseent where it doesn't 

-Requires the Hessian 
-Needs to solve a linear system of 
equations 

Modified Newton 

-Maintains a deseent property 
-Avoids divergenee 

-Requires a line seareh 


Table 1. Methods of Multi-Variable Uneonstrained Optimization 


3. Constrained Optimization 

Constrained optimization is the minimization of an objeetive funetion subjeet to 
the eonstraints of the possible values of the independent variables. Constraints ean be 
linear or nonlinear with either equality or inequality eonstraints. 

The typieal equality eonstrained optimization problem has the form: 

minimize f(x) 
subjeet to h(x) = 0 

where f(x) is the sealar valued funetion and h(x) is the veetor valued eonstraint 
funetion. Beeause of the eonstraint, stationary points (minimum, maximum, or saddle 
points) of f(x) may not be valid solutions sinee they may not satisfy the eonstraints. The 
eommon approaeh to solving eonstrained optimization problems is with the method of 
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Lagrange multipliers. This method converts the constrained optimization problem into an 
unconstrained one, which is then solved with one of the methods listed on Table 1. The 
stationary points of the Lagrangian are potential solutions of the constrained optimization 
problem [Ref 12], 

J=i 

The typical inequality constrained optimization problem has the form; 

minimize f(x) 
subject to g(x) < 0 

where f(x) is the scalar valued function and g(x) is the vector valued constraint 
function. The discussion of inequality constraints is lengthy and thoroughly discussed in 
[Ref 11] and [Ref 12]. 

E, DECISION SUPPORT SYSTEMS 

Decision Support Systems (DSS) are a specific class of computerized information 
systems that support the decision makers. DSS are typically interactive computer-based 
systems intended to aide decision makers by using communications technologies, data, 
documents, knowledge and/or models to identify, solve problems, and make decisions. 
The type of DSS developed for this thesis is a model-driven DSS. The model-driven 
DSS emphasizes access and manipulation of optimization and simulation models. It uses 
data and parameters provided by the designer, are not usually data intensive, and do not 
use large databases. 


20 



III. INVESTIGATION AND ENHANCEMENTS 


A. CORRECTIONS TO EXISTING NUMERICAL ISOLATOR MODELS 

The Bilinear, Maxwell, and Wen models were adopted from [Ref. 8]. In this 
referenee, the author presents their use with TSS and RBB to solve the base isolation and 
struetural response for three struetural models. The MATLAB programs inherited from 
this researeh failed to duplieate the results presented in [Ref. 8]. The isolator model 
programs were thoroughly investigated for model theory, aeeuraey, program flow, and 
eompatibility with RBB. The results of this investigation are presented in Table 2. Other 
errors of minor importanee inelude: 

• The El Centro time history plot on [Ref. 8, Page 53] is ineorreetly presented in 
inehes. The maximum and minimum displaeement points are eonsistent with the 
measurements in eentimeters. The data presented in this referenee eontains unit 
dimension errors. 

• The starting veetor of isolator forees used for TSS was a veetor of zeros. A better 
veetor of starting values is a veetor of ones, as reeommended in [Ref. 14]. 

If a simulation is performed with only one bloek, RBB is not used. Any isolator 
model not written to “remember” forees and displaeements from previous bloeks will 
funetion eorreetly. Onee the number of bloeks exeeeded one, RBB is employed and 
isolator models that don’t “remember” will report a zero foree value at the beginning of 
eaeh bloek. Examples of these foree diseontinuities are show in Eigures 15 - 17. The 
program-logie flowehart for proeessing RBB bloeks is presented in Eigure 18. It is 
imperative that MATEAB’s persistent variables, those that retain the final foree and 
displaeement values for the previously proeessed bloek, be elear from memory. 
Otherwise the persistent variables will remain in memory until the user terminates the 
MATEAB session. Any further proeessing without first elearing these variables will 
result in ineorreet results. 
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Isolator Model 

Problems 

Corrections 

Bilinear 

• No funetionality with RBB. 

• Foree diseontinuities for more 
than one block. 

• Poor program flow. 

• Poorly documented. 

• Model re-written to function 
properly with RBB. 

• Logic and program flow 
improved - properly documented 

• Wrote a separate model that does 
not account for strain hardening 
since LRB are self annealing at 
standard temperature. 

Maxwell 

• No functionality with RBB. 

• Force discontinuities for more 
than one block. 

• Runge-Kutta 4* order method 
solved incorrectly for range of 
time steps. 

• Poor program flow. 

• Poorly documented. 

• Model re-written to function 
properly with RBB. 

• Logic and program flow 
improved 

• Runge-Kutta 4* order method 
solves correctly for any user- 
defined time step. 

• Logic and program flow 
improved - properly documented 

Wen 

• Incorrect Wen equation 

• No functionality with RBB. 

• Force discontinuities for more 
than one block. 

• Runge-Kutta 4* order method 
solved incorrectly for range of 
time steps. 

• Poor program flow. 

• Poorly documented. 

• Model re-written with correct 
equations [Ref 13] and to 
function properly with RBB. 

• Logic and program flow 
improved 

• Runge-Kutta 4* order method 
solves correctly for any user- 
defined time step. 

• Logic and program flow 
improved - properly documented 


Table 2. Problems Diseovered with Numerieal Isolator Models 
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Figure 15. Force Discontinuities in with the Bilinear Model (Force vs. Disp) 



Figure 16. Force Discontinuities with the Maxwell Model (Force vs. Time) 
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Figure 17. Force Discontinuities with the Maxwell Model (Force vs. Disp) 



Figure 18. Blocked Data Flow for an Isolator Model that Requires Retained Values 
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1, Wen Model Validation 

The Wen isolator model is of partieular interest due to its ability to represent all 
varieties of hysteretie behavior by the proper seleetion of parameters for its governing 
equation. [Ref. 13] The restoring force for this hysteretie damper is; 

F, 

Q{x,x) = a—x + {\-a)F^,Z (10) 


Here, Z is the dimensionless hysteretie component, which satisfies the following 
nonlinear first order differential equation; 

YZ = -y|ji;|z|zr‘-y0ji;|zr+Hji: (11) 


The Wen numerical model program developed for this thesis accurately 
reproduced the results presented by Constantinou and Tadjbakhsh [Ref 13]. The 
following figure shows the duplicated results for the base displacement of 5.5 inches with 
a frequency of 0.9 Hz used for their study. 




Figure 19. Hysteresis Loops for Different Wen Parameters 
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Figure 20. Hysteresis of: (a) Torsional Device; (b) LRB from Ref [13] 

For a detailed review of the dimensionless parameters a, P, y, A, and r\, please 
refer to [Ref. 13]. 

B, SLIDING NUMERICAL ISOLATOR MODELS 

Sliding bearings will limit the transmission of lateral forces to the isolated 
structure. The three sliding models presented here were adopted from [Ref 7]. This 
reference is maintained by the Federal Emergency Management Agency (FEMA) and 
sets guidelines the analysis and design of passive seismic isolation and energy dissipation 
systems. These models do not take into account the effects of slipping and sticking, 
which can transform low-frequency base vibrations into high-frequency structural 
vibrations [Ref 2]. 

1. Flat Sliding Surface with Restoring Device 

This friction model resembles a linear spring in parallel with a linear damper. 
The governing equation for this model is stated in Equation (10) and does not require 
retained values from previous blocks to function correctly with RBB. 

F = kx + jU^N sgn{x) (12) 
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This equation uses the eoeffleient of statie frietion, qs, sinee the isolator is at rest 
before a disturbanee. Also, the dynamie frietion eoeffieients for materials eommonly 
used with sliding bearings differs very little from the statie frietion eoeffieient [Ref. 2], 
This model aeeepts the base displaeement, base veloeity, normal foree (the distributed 
foree of the strueture on eaeh isolator), and a frietion eoeffieient in order to proeess and 
return the isolator restoring foree. Flat sliding bearings with restoring eapabilities are 
usually the eombination of a sliding and elastomerie bearing. 

2, Spherical Sliding Surface 

In this model the spherieal sliding surfaee aets as the restoring deviee. Notiee that 
as the sphere’s radius of eurvature, R, inereases, the restoring foree part of Equation (11) 
beeomes irrelevant and degrades to a flat sliding bearing with out a restoring foree. This 
model does not require retained values from previous bloeks to funetion eorreetly with 
RBB. 


F = 


A 

R 


x + ju^Nsgn(x) 


(13) 


3, Conical Sliding Surface 

The eonieal sliding surfaee model does require the last foree value from the 
previous bloek in order to funetion eorreetly with RBB. This is due to the transition 
between the eone’s spherieal surfaee and the inelined surfaee. Note that the foree vs. 
displaeement plot from Figure 9 shows two possible foree values for the same 
displaeement. Retention of the previous foree is needed to determine if the displaeement 
is eausing an inerease in foree (top line) or a deerease (bottom line). The program uses 
Equation (14), whieh was derived from plane geometry, to determine the transition 
between the spherieal surfaee and the inelined surfaee. 


s 


transition 


1 


1 

ian{(/)Y 


+ 1 


(14) 


If the eriteria of -Stmnsition <x< Stransuion is met, then the program uses Equation 
(14) ealeulate the resulting foree. Outside of the transition boundaries, the program uses 
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Equation (15) to calculate the restoring force. Note that the value of Stransiuon decreases as 
the incline angle decreases, which makes Equation (13) the dominating one. 

^ = ^„±2//,iVtan(^) (15) 



Eigure 21. Conical Sliding Surface Numerical Model (Eorce vs. Disp) 

C. CONVOLUTION-DRIVEN NUMERICAL INTEGRATION 

The programs that run TSS and RBB exploit MATLAB’s conv command for 
numerical integration. The convolution theorem states that the convolution of two 
vectors is the same as multiplying their Eourier transforms. MATLAB's Signal 
Processing Toolbox is optimized for speed and can perform convolutions very efficiently. 
The following convolution example should help the reader understand and appreciate 
convolution-driven numerical integration. 

To integrate the liney = x, let x:= [0 1 2 3 4 5] and v = [1 1 1 1 1 1]^. Before 
convolving these two vectors, we apply a filter matrix to x and discard the last four 
columns. 
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( 16 ) 


(17) 


The answer is = [0 1 3 6 10 15]^, whieh is the reetangular rule integration of 
y=x. The analytieal answer is jo = [0 1 4 9 16 25]^. In order to improve the eonvolution 
aeeuracy, we weight v with the Composite Trapezoid Rule, therefore v = [1 2 2 2 2 1]^. 
After a matrix multiplieation as in Equation (15), the new solution is j/ = [0 1 4 9 16 25]^, 
whieh for this ease happens to be the exaet answer. This method allows the easy 
applieation of any integration rule by using a veetor of weightings eharaeteristie to the 
integration rule desired. The only requirement for this method to work eorreetly is that 
the square matrix must be lower triangular and the diagonal must be zero. This is not a 
problem with TSS sinee the diagonal of the Impulse Response Funetion (IRE) matrix is 
always zero at t = 0. The following is a list of integration rules with their respeetive 
eharaeteristie weighting veetor. 

1. Trapezoid Rule - Error Oh^ 

J /(x) Jx = [/(Xq ) + 2/(Xi) + /(X3)] (18) 


v = -(l 2 2 ... 2 2 1) 


(19) 
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( 20 ) 


2, Simpson’s Rule - Error Oh^ 

I /(x) Jx = y [/(Xq ) + 4/(Xi) + 2/(X2 ) + 4/(X3 ) + /(X 4 )] 

v = ^(l 4 2 4 2 ... 2 4 2 4 if (21) 

3, Weddle’s Rule - Error Oh^ 

h 

J f{x)dx = -[41/(x„) + 216/(x,) + 27/(xJ + 212f{x^) + 27/(xJ + 216/(x,) + 41/(xJ] (22) 

a 140 

v = —(41 216 27 272 27 216 ... 216 27 272 27 216 41)^ (23) 

140^ 7 V / 

When using these integration weightings with TSS, the length of v must be the 
same length as the veetor of isolator synthesized forees. The following example will 
demonstrate how RBB works and how to ineorporate the veetor of integration 
weightings. Figure 22 presents a proeess in its fourth iteration with the following 
eonditions; 

• This data was divided into four equally-sized blocks and manipulated as in 
Equation (17). 

• The initial guess was/= [11111111111]^. 

• First three blocks have converged and/was updated. 
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Figure 22. Example of a Recursive Block-by-Block Process 


Each successive block will solve faster since it relies on the synthesized answers 
from prior blocks. Before the fourth block can be synthesized, the following pre-process 
must take place. 
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v6y 


• update x(10; 12) 


• update x(10:12) 


update x(10; 12) 


(24) 


The fourth /block and the updated fourth x block are then processed by the TSS 
algorithm. To improve the accuracy of RBB, apply the blocked vector of integration 
weightings, IW, as shown in Eigure 23 - pad with zeros to match the length of f 
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Figure 23. Example of a Reeursive Block-by-Block Process with Simpson’s Rule 


The vector of integration weightings is element-wise multiplied with the vector of 
forces. Note that the padded zero has no effect on the answer since the matrix diagonal is 
already zero. Figure 24 shows the accuracy improvements for a 2-DOF spring-mass- 
damper system solved with TSS / RBB (8 blocks) and Simpson’s Rule integration. The 
results were compared to a Working Model 2-D simulation with the same parameters. 



Figure 24. Accuracy Improvement for a 2-DOF SMD System 
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Deviation from Working Model Results 



Figure 25. Differential Integration Error for 2-DOF SMD System 

D. PROGRAM OPTIMIZATION FOR SPEED 

The MATLAB programs presented here were thoroughly investigated, modified, 
and/or re-written to remove proeess redundaneies and improve program flow and inter- 
eompatibility. Fleavy use was made of MATLAB’s veetorization and eonvolution 
eapabilities, whieh are optimized for speed and eommonly used for signal proeessing. 
The Wen and Maxwell isolator models employ a Runge-Kutta 4**' order method to solve 
the governing differential equations (GDE). The GDEs eould have been solved faster 
with MATLAB’s ODE 45 solver, but it was deeided that the present method, whieh uses 
EOR; NEXT loops, would suffiee. The eompiled version of this program, like other 
programs, will exeeute faster than in its native programming environment. 
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E. OASIS (GUI) 

Optimization and Analysis of Structural Isolation Systems (OASIS) is the DSS 
tool developed to bring together and automate the optimization and analysis process. The 
following figures show screen shots of the GUI modules. 



Figure 26. Main GUI Module 
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Figure 27. Transient Structural Synthesis Module 



Figure 28. Isolator Options and Constraints 
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B dit OptimizMion Paran^ter 

—Model Selected for Optimization- 


7DQF1.syn 


Select Input File 


Select Output Data File 


—Select Isolation Type- 

[Bilinear w/oStrainHardenin 

Edit Constraints | 

—Minimize- 

Analysis Point Displacement 

Analysis Point Acceleration 


f—Maximum- 


I solator D isplacement (R elalive) ^ 
Analysis Point Acceleration T 


-Coupling Point(s)- 



0.52 mm 

mm/'s''2 


1000 


I— Analysis Point(s) - 



Edit Optimizer Options 
Start Optimization 


B FMINCON Options 


—FMINCON Optimset Parameters- 

DerivativeCheck F* 

Diagnostics 

r 

DiffMaxChange 

1 o-i 

DiffMinChange 

1 1e-008 

Display 

liter 

G oalsE xactAchieve 

1 0 

GradConstr 

r 

GradObj 

r 

Hessian 

r 

HessMull 

1 - 

HessUpdate 

|bfgs 

Jacobian 

r 

JacobMult 

1 - 

LargeScale 

r 

LineSearchType 

] quadcubic 

LevenbergMarquardt 

r 

MaxFunEvals 

1 500 

Maxi ter 

1 30 

MeritFunction 

1 singleobj 

PrecondBandWidlh 

1 0 

T olCon 

1 le-OO? 

T olFun 

1 le-OO? 


Save and Close Restore Defaults 


Figure 29. Optimization Module with Advanced FMINCON Options 
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IV. 


PROGRAM VERIFICATION 


A, WORKING MODEL 2-D 

Working Model 2-D is a conceptual design tool that allows users to create and 
analyze real-life mechanical systems. It was extensively used to verify the accuracy of 
the results generated with OASIS. Presented in the following tables and figures are the 
validation results for 2-DOF and 7-DOF models with parallel spring-dampers 
respectively. Cset is a vector with the index value of the coupling node(s). For example, 
the following 2DOF system is coupled at mi (node #1), therefore Cset = [1]. 


1, 2DOF System 

Problem: Determine the response of m 2 . 

Model: mi = 1 slug 

m 2 = .5 slug 
k = 50 Ibf / ft 
c =0 Ibf-s / ft 


Isolator: 

Spring - Damper in parallel 
k= 10 Ibf/ft 
c = 0.5 Ibf-s / ft 

Cset: 

[1] 

Base Motion: 

Sinusoidal 

A = 0.0833 in 
/= 10 Hz 

Duration: 

2 sec. 

Results: 

Baseline results from WM2D. 

Max Disp: 

0.0037 ft 

Statistics: 

Synthesis deviation from WM2D. 

Average: 

9.6100 X lO'"’ 

Standard Deviation: 2.7574 x 10'^ 

Maximum Error: 6.3873 x 10'^ 
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2DOF - Response for Mass #2 



Figure 30. 2DOF Validation Model - Response for Mass #2 


2, 7DOF System 

Problem; Determine the response of m 2 . 

Model: mi = 3.5 slug 

m 2 = 3 slug 
m 3 = 2.5 slug 
m 4 = 2 slug 
ms = 1.5 slug 
me = 1 slug 
m? = 0.5 slug 

ki = k 2 = k 3 = k 4 = k 5 = k 6 = 50 Ibf / ft 
ei = e 2 = C3= C4= C5= C6= 0 Ibf-s / ft 


Isolator; Spring - Damper in parallel: 

k = 200 Ibf/ft 
c = 2.5 Ibf-s / ft 

Cset; [1 7] 

Base Motion; Sinusoidal 

A = 0.08333 ft 
/= 10 Hz 

Duration; 2 see 


Results; Baseline results from WM2D. 

Max Disp; 0.0054 ft 
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Statistics; Synthesis deviation from WM2D. 
Average; 8.9500x 10'^ 

Standard Deviation; 3.6670 x 10'"^ 

Maximum Error; 9.6700 x 10'^ 


7DOF - Response for Mass #2 



Figure 31. 7DOF Validation Model - Response for Mass #2 


B. SINGLE BAY - FOUR STORY BUILDING FRAME 

Presented here is a validation of a four story building isolated with Maxwell 

elements. The force vs. displacement plot of the Maxwell isolator, [Ref. 8, Page 59], is 

presented correctly, but was most likely processed in one block since no force 

discontinuities are present. The following validation was preformed to reproduce the 

results from [Ref. 8] with the RBB algorithm. 

1, Four Story Building Frame 

Problem: Determine a corner lateral displacement 

Model; bldgdstoryfreehor.mat 

Isolator; Maxwell Isolator 

k = 300 Ibf / in 
c = 51 Ibf-s / in 
Cset; [1 2 3 4] 

Base Motion; ElCentrol940NS.mat 
Duration; 50 sec. 
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Figure 32. Lateral Flysteresis of Corner Node , 16 bloeks 



Isolator 

Figure 33. Flysteresis of Corner Node (From Ref. 8) 
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V. SAMPLE OPTIMIZATION PROBLEMS 


A. SINGLE BAY FOUR-STORY BULBING 

The four-story building used for [Ref 8] is revisited here with a variety of 
isolators. This strueture is assembled with 50 ksi steel struetural members with the 


following specifieations; 

• Columns: W36x486 

• First floor beams: W36xl70 

• Second floor beams: W36xl70 

• Third floor beams: W36xl60 

• Fourth floor beams: W36xl50 

• Roofbeams: W36xl35 


The following objective functions use the absolute values of the design variables 
to ensure the minimum values are not interpreted as large negative values. 

1, Parallel Spring and Damper 


Base Motion 

Sine, lin, 10 Hz, 7 sec 

Coupling Nodes (cset) 

[1 2 3 4] 

Analysis Node (rset) 

[20] - 5*'' floor node 

Starting Isolator Values 

k = 15 Ib/in ; c = 0.5 Ib-s/in 

Objective 

Minimize 5**' floor lateral displacement (node 15) 

Optimal Isolator Values 

k = 15 Ib/in ; c = 0.05 Ib-s/in 

Max Displacement 

Max 5* floor lateral disp = 0.0032 in 

Process Time 

00:09:42 


Table 3. Optimization Results for 4-Story Building with Parallel Spring - Damper 

Isolators 
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Objective Function; 


MIN; 

X20 


(5*’’ floor lateral displacement) 

Subjeet to; 

k-30 

<0 

(Isolator constraints) 


-k+ 10 

<0 



c-0.7 

<0 



-e + O.Ol 

<0 



-X20 

<0 

(Nonnegativity) 



Figure 34. Parallel Spring - Damper; Foree vs. Displaeement - Starting Values 



Figure 35. Parallel Spring - Damper; Foree vs. Displaeement - Optimal Values 
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)( 1 Q'^ Response Motion for element ffiO 



Figure 36. Parallel Spring - Damper; Displaeement of node #20- Before Optimization 


x iq'^ Response Motion for element #20 



Figure 37. Parallel Spring - Damper; Displacement of node #20- After Optimization 
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2. 


Flat Slider with Restoring Spring 


Base Motion 

Sine, lin, 10 Hz, 7 sec 

Coupling Nodes (eset) 

[1 2 3 4] 

Analysis Node (rset) 

[20] - 5*’’ floor node 

Starting Isolator Values 

k = 10 Ib/in ; p = 0.5 ; N = 2183 lb 

Objeetive 

Minimize 5**' floor lateral acceleration 

Optimal Isolator Values 

k = 23.3 Ib/in ; p = 0.1 

Max Displaeement 

Max 5**^ floor lateral accel = 13.0446 in/sec^ 

Proeess Time 

00:03:50 


Table 4. Optimization Results for 4-Story Building with Flat-Slider Restoring Spring 

Isolators 

Objective Function; 


MIN: 

Vo 


(5**' floor lateral acceleration) 

Subject to: 

k-30 

<0 

(Isolator constraints) 


-k + 5 

<0 



o 

1 

<0 



-p + O.I 

<0 



- XjQ 

<0 

(Nonnegativity) 
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Figure 38. Flat Slider with Restoring Spring; - Before Optimization 
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Figure 39. Flat Slider with Restoring Spring; - After Optimization 


47 




















































Figure 40. Flat Slider with Restoring Spring: Lateral Acceleration of node #20 - Before 

Optimization 



Figure 41. Flat Slider with Restoring Spring: Lateral Displacement of node #20 - After 

Optimization 
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3, Spherical Slider 


Base Motion 

Sine, lin, 10 Hz, 7 see 

Coupling Nodes (eset) 

[1 2 3 4] 

Analysis Node (rset) 

[20] - 5*’’ floor node 

Starting Isolator Values 

p = 0.5 ; r = 60 in; N = 2183 lb 

Objeetive 

Minimize 5**' floor lateral displaeement 

Optimal Isolator Values 

p = 0.1; r = 120 in 

Max Displaeement 

Max 5**^ floor lateral disp = 0.0021 in 

Proeess Time 

00:04:06 


Table 5. Optimization Results for 4-Story Building with Spherieal-Slider Isolators 

Objective Function; 


MIN; 

X20 


(5**' floor lateral displaeement) 

Subjeet to: 

p - 30 

<0 

(Isolator eonstraints) 


-p + 5 

<0 



r-0.7 

<0 



-r + O.I 

<0 



-X20 

<0 

(Nonnegativity) 
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Figure 42. Spherical Slider; Force vs. Displacement - Before Optimization 



Figure 43. Spherical Slider; Force vs. Displacement - After Optimization 
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Figure 44. Spherical Slider; Lateral Displacement of node #20- Before Optimization 



Figure 45. Spherical Slider; Lateral Displacement of node #20- After Optimization 
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4, Bilinear - No Strain Hardening 


Base Motion 

Sine, lin, 10 Hz, 10 sec 

Coupling Nodes (eset) 

[1 2 3 4] 

Analysis Node (rset) 

[20]- 5*^ floor node 

Starting Isolator Values 

ki= 4937 Ib/in; yp = 0.4; Fyt = 1997 lb 

Fye = -1993 lb 

Objeetive 

Minimize 5*' floor lateral displacement 

Maximum relative isolator displacement = 2 in 

Optimal Isolator Values 

ki= 4933 Ib/in; yp = 0.25; Fyt = 2100 lb 

Fye = -1900 1b 

Max Displacement 

Max 5*' floor lateral disp = 0.041 in 

Process Time 

00:07:12 

Bearing Data 

High-Damping Natural Rubber 

Dia: 19.68 in 

Length: 19.68 in 

Shear Mod: 79.77 Psi 

Stiffness: 4933.21 lb / in 


Table 6. Optimization Results for 4-Story Building with Bilinear NSH Isolators 
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Objective Function; 


MIN: 

X 20 


(5**' floor lateral displacement) 

Subject to; 

ki -5200 

<0 

(Isolator constraints) 


-ki+ 4800 

<0 



yp - 0.5 

<0 



-yp + 0.1 

<0 



Fyt-2200 <0 

-Fyt + 1800 < 0 
Fyc - 2200 < 0 

-Fyc + 1800 < 0 

-X 20 < 0 (Nonnegativity) 
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Figure 46. Bilinear NSFI; Flysteresis Plot - Before Optimization 



Figure 47. Bilinear NSH; Flysteresis Plot - After Optimization 
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Response Motion for element #20 



Figure 48. Bilinear NSFI; Lateral Displacement of node #20 - Before Optimization 



Figure 49. Bilinear NSFI; Lateral Displacement of node #20 - After Optimization 
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5. 


Maxwell 


Base Motion 

Sine, 12in, 5 Hz, 10 sec 

Coupling Nodes (cset) 

[1 2 3 4] 

Analysis Node (rset) 

[1 2 3 4] - Isolator Node 

Starting Isolator Values 

k = 3000 Ib/in; c = 50 Ib-s/in 

Objective 

Minimize isolator node absolute displacement 

Optimal Isolator Values 

k = 3000 Ib/in; c = 40 Ib-s/in 

Max Displacement 

Max isolator ode absolute disp = 0.08 in 

Process Time 

00:01:55 


Table 7. Optimization Results for 4-Story Building with Maxwell Isolators 


Objective Function; 

MIN: 
Subject to: 


Xi + X 2 + X 3 + X 4 

k-3000 =0 

c - 60 < 0 

-c + 40 < 0 

Xi, X2, X3, X4 > 0 


(Isolator lateral displacements) 
(Isolator constraints) 

(Nonnegativity) 
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){'IQ*' Force vs. Relative Displacement for isolator # 1 



Figure 50. Maxwell; Force vs. Displacement Plot - Before Optimization 



Figure 51. Maxwell; Force vs. Displacement Plot - After Optimization 
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Figure 52. Maxwell; Lateral Displaeement - Before Optimization 



Figure 53. Maxwell: Lateral Displacement - After Optimization 
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6, Wen - Sine Base Motion 


Base Motion 

Sine, 12m, 5 Hz, 10 see 

Coupling Nodes (eset) 

[1 2 3 4] 

Analysis Node (rset) 

[1 2 3 4] - Isolator Nodes 

Starting Isolator Values 

Yield Foree = 3000 lb 

Yield Displaeement(Y) = 3 in 

Yield to Post (a) = 0.1 

P = 5 

r = 5 

A = 3 

r| = 3 

Objeetive 

Minimize isolator node absolute displaeement 

Optimal Isolator Values 

Yield Foree = 2970 lb 

Yield Displaeement(Y) = 10.8 in 

Yield to Post (a) = 0.1 

p= 100 

r = 313.4 

A = 2 

r| = 3 

Max Displaeement 

Max isolator node absolute disp = 0.1350 in 

Proeess Time 

02:01:53 


Table 8. Optimization Results for 4-Story Building with Wen Isolators - Sine 

Displaeement 
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Objective Function; 


MIN; 

Xi + X2 + X3 + X4 

(Isolator lateral displacements) 

Subject to; 

Fy -3500 

<0 

(Isolator constraints) 


- Fy + 2500 

<0 



Y- 100 

<0 



- Y + 2.7 

<0 



a - 0.5 

<0 



- a + 0.05 

<0 



p- 1000 

<0 



-p + 3 

<0 



F- 1000 

<0 



-r + 3.5 

<0 



A-5 

<0 



-A+ 1 

<0 



0-5 

<0 



-r| + 1 

<0 



Xi, X2, X3, X4 

>0 

(Nonnegativity) 


ri, A 

=Integer 
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Force vs. Relative Displacement for isolator#! 



Figure 54. Wen; Flysteresis Plot - Before Optimization 


Force vs. Relative Displacement for isolator # 1 



Figure 55. Wen; Flysteresis Plot - After Optimization 
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Response Motion for element #1 



Figure 56. Wen; Lateral Displaeement of node #1- Before Optimization 


Response Motion for element # 1 



Figure 57. Wen; Lateral Displacement of node #1 - After Optimization 
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7. Wen - Kobe Japan 1995 NS 


Base Motion 

KobeJapanl995NS US.mat20sec 

Coupling Nodes (eset) 

[1 2 3 4] 

Analysis Node (rset) 

[20] - 5*’’ floor node 

Starting Isolator Values 

Yield Force = 3000 lb 

Yield Displacement(Y) = 3 in 

Yield to Post (a) = 0.1 

P = 5 

r = 5 

A = 3 

r| = 3 

Objeetive 

Minimize 5**^ floor lateral displacement. 

Isolator displacement relative to ground motion 

cannot exceed 10 in. 

Optimal Isolator Values 

Yield Force = 2968 lb 

Yield Displacement(Y) = 5.9 in 

Yield to Post (a) = 0.5 

p = 98.1 

r = 313.4 

A = 2 

r| = 3 

Max Displaeement 

Max 5**^ floor node lateral disp. = 3.193 in 

Max isolator node relative disp = 9.1371 in 

Process Time 

02:05:44 


Table 9. Optimization Results for 4-Story Building with Wen Isolators - 

KobeJapanl995NS_US.mat 
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Objective Function; 

MIN: X 20 

Subject to: Fy -3500 

- Fy + 2500 
Y- 100 

- Y + 2.7 
a - 0.5 

- a + 0.05 
p- 1000 
-p + 3 
F- 1000 
-r + 3.5 
A-5 
-A+ 1 
0-5 

-r| + 1 

Zi^ Z2, Z3^ Z4 

X20 

ri, A 


(5* floor lateral displacements) 

< 0 (Isolator constraints) 

<0 

<0 

<0 

<0 

<0 

<0 

<0 

<0 

<0 

<0 

<0 

<0 

<0 

<10 in (Isolator relative displacements) 

> 0 (Nonnegativity) 

= Integer 
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Figure 58. Wen; Flysteresis Plot - Before Optimization 



Figure 59. Wen; Flysteresis Plot - After Optimization 
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Relative Response Motion for Element #1 



Figure 60. Wen; Relative Lateral Displacement of node #1 - Before Optimization 


Relative Response Motion for Element #1 



Figure 61. Wen; Relative Lateral Displacement of node #1- After Optimization 
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Response Motion for element #20 



Figure 62. Wen; Lateral Displaeement of node #20 - Before Optimization 



Figure 63. Wen; Lateral Displaeement of node #20- After Optimization 
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Figure 64. Wen; Lateral Aeeeleration - Before Optimization 



Figure 65. Wen Lateral Aeeeleration - After Optimization 
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Figure 66. No Isolation; Lateral Displacement of node #20 



Figure 67. No Isolation; Lateral Acceleration of node #20 
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B, FIFTEEN BAY THIRTY-STORY BULBING 

The thirty-story building used for [Ref 8] is revisited here with Wen isolation. 
This strueture is coupled with 81 isolators and exposed to both sinusoidal and Kobe 
earthquake base motions. This model was processed to demonstrate OASIS's ability to 
quickly process large-complex structural systems. 


1. Wen - Sine Base Motion 


Base Motion 

Sine, 12in, 5 Hz, 10 sec 

Coupling Nodes (cset) 

[1 2 3 4 ... 77 78 79 80 81] (81 nodes) 

Analysis Node (rset) 

[1215 2511]- 15'%nd 3 floors 

Isolator Values 

Yield Force = 2970 lb 

Yield Displacement(Y) = 10.8 in 

Yield to Post (a) = 0.1 

p= 100 

r = 313.4 

A = 2 

r| = 3 

Objective 

Determine the response of the 15* and 3D* floor 

nodes. 

Max Displacement 

Max Lateral Disp 15* Floor; 0.0246 in 

Max Lateral Disp; 3D* Floor: 0.0251 in 

Process Time 

03:37:43 


Table 10. Synthesis Results for30-Story Building with Wen Isolators - Sine Displacement 
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Figure 68. Thirty-Story Building Wire Frame (From Ref. 8) 


72 





Figure 69. 30-Story Building; Isolator #1 Lateral Displaeement 


Force vs. Relative Displacement for Isolator#! 



Figure 70. 30-Story Building; Isolator #1 Hysteresis Plot 
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Response Motion for IS*'^ Floor Node 



Figure 71. 30-Story Building; 15*^ Floor Lateral Displacement 



Figure 72. 30-Story Building; 3 L* Floor Lateral Displacement 
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Figure 73. 30-Story Building; 31'^* Floor Acceleration 
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2, Wen - Kobe Japan 1995 US Base Motion 


Base Motion 

KobeJapanl995NS US.mat20sec 

Coupling Nodes (eset) 

[1 2 3 4 ... 77 78 79 80 81] (81 nodes) 

Analysis Node (rset) 

[1215 2511]- 15*%nd 3 floors 

Isolator Values 

Yield Force = 2970 lb 

Yield Displacement(Y) = 10.8 in 

Yield to Post (a) = 0.1 

p= 100 

r = 313.4 

A = 2 

r| = 3 

Objeetive 

Determine the response of the 15* and 3D* floor 

nodes. 

Max Displaeement 

Max Lateral Disp 15* Floor; 3.6817 in 

Max Lateral Disp: 3D* Floor: 3.6841 in 

Process Time 

06:27:36 


Table 11. Synthesis Results for 30-Story Building with Wen Isolators - 

KobeJapanl995NS_US.mat 
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Figure 74. 30-Story Building; Isolator #1 Lateral Displacement 


Force vs. Relative Displacement for Isolator #1 



Figure 75. 30-Story Building; Isolator #1 Hysteresis Plot 
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Figure 76. 30-Story Building; 15* Floor Lateral Displaeement 



Figure 77. 30-Story Building: 3 Floor Lateral Displacement 
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Figure 78. 30-Story Building; 3Floor Acceleration 

2 

(Time (sec) vs. Acceleration (in/sec ) 
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Time vs. Lateral Displacement of a 31®* Floor Node 



Figure 79. No Isolation; 31*^ Floor Lateral Displaeement 


Time vs. Lateral Acceleration of a 31®* Floor Node 



Figure 80. No Isolation; 3 L* Floor Lateral Acceleration 
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VI. CONCLUSION 


MATLAB programs for Transient Structural Synthesis (TSS), Recursive Block- 
by-Block (RBB), and nonlinear numerical isolator models from [Ref 8] were thoroughly 
investigated. As a result, numerous discrepancies in the isolator models were corrected, 
program logic was optimized, and a method that transforms MATLAB’s rectangular-rule 
convolution into any integration rule with higher-order error was developed. This 
research also led to the development of three numerical models for sliding-friction 
isolators and a user-friendly Decision Support System (DSS) tool for transient response 
and analysis. OASIS (Optimization and Analysis of Structural Isolation Systems) is a 
DSS tool that gives a user the ability to quickly change isolator parameters for TSS / 
RBB, change types of base motions, optimize the isolator parameters, and post-synthesis 
/ optimization analysis tools. 

The program results were carefully validated with 2-DOF and 7-DOF systems in 
Working Model 2-D (WM2D). The results from OASIS and Working Model 2-D 
compared very well for short simulation times. For longer simulation times (t > 15 sec), 
small deviations in response were noted. This is due to Working Model's attempt to 
solve a "real world" system of blocks, springs, and dampers, instead of a lumped mass 
system. OASIS uses mass-normalized mode shapes in its calculation, which result in 
better solutions when compared to WM2D. 

Several optimization scenarios were presented for comparison. Each example 
used different isolators with initial, lower-bound, and upper-bound parameters. The user 
also has the option to "fix" certain parameters with the use of nonlinear equality 
constraints. The objective function is the linear sum of either displacements, 
accelerations, or both. However, a future version will allow the optimization of nonlinear 
objective functions. 

OASIS can also be used to analyze and optimize models of smaller scale. The 
program only needs the mode shapes and natural frequencies for the coupling nodes and 
the analysis nodes. TSS makes it possible to reduce a large DOF system into a small one 
that can be solved quickly. For example, the four story structure used in the optimization 
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examples has 20 nodes, eaeh with six degrees of freedom. Only five nodes with lateral 
degrees of freedom are needed to find the lateral response of a node on the fifth floor - 
four eoupling nodes and one fifth floor node. OASIS's ability to aeeept any size model, 
eouple it with an isolator, and apply any base motion greatly extends its use for optimal 
engineering design. 
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VII. RECOMMENDATIONS 


OASIS does not account for the uplift effects created by a structure's rocking 
during a base-motion event. To consider uplift in the structure's response, the vertical 
degrees of freedom need to be included for the impulse response function formulation. 
Reference [8] presents results of programs developed to account for uplift by modeling a 
linear spring. More information on past and current research on the effects of uplift and 
base isolation can be found in [Ref. 17]. 
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